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Dynamics of fluctuations in unstable Bose-Einstein con- 
densates is analyzed by the solution of approximate operator 
equations. In the case of a condensate with a negative scatter- 
ing length the present treatment describes a delay of collapse, 
in agreement with recent experiments. In the case of a colli- 
sion of two condensate wavepackets it is shown that quantum 
effects lead to a Bose enhancement of elastic-scattering losses. 
In both cases the noncondensate atoms are formed as entan- 
gled pairs in squeezed states. 
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I. INTRODUCTION 

Recent progress in the physics of Bose-Einstein con- 
densation (BEC) of atomic gases (see Refs. and ref- 
erences therein) stimulated a further development of the- 
oretical methods used previously for the description of 
superfluidity, superconductivity, and collective effects in 
heavy atoms and nuclei (see Refs. These methods 

treat condensates as mean fields, while interactions be- 
tween particles, according to the Bogolubov theory, lead 
to formation of fluctuations (see Refs. @-|l^]). In the 
case of a positive scattering length, the condensate (mean 
field) is stable, while fluctuations can be treated by using 
a Bogolubov transformation to quasiparticle states. 

Recent experiments |Tl|-|l4| drew attention to the 
physics of condensates with a negative scattering length. 
In this case, a free condensate (with a homogeneous mean 
field) is unstable, collapsing within a finite time ||, while 
a trapped condensate is stable at the low occupations of 
the experiment |OJI3|. In the JILA experiments on Rb 
JL^JTijl , the effect of Feshbach resonance (see Ref. 15 



was used for changing the sign of the scattering length. 
This effect allows the formation of a condensate with a 
positive scattering length, followed by a study of the dy- 
namics of collapse after applying a sudden change of sign 
of the scattering length. A formal extension of the Bo- 
golubov transformation to the case of negative scattering 
lengths would lead to complex energies of the quasipar- 
ticles, but actually such a transformation does not exist 
in this case, as was proven by Bogolubov (see Ref. ||). 
Nevertheless, the development of fluctuations can sub- 
stantially affect the dynamics of an unstable condensate 
(see the theoretical analysis [[l6| of the JILA experiments 



The present paper describes the dynamics of fluctua- 
tions by using the parametric approximation employed 
in Ref. |17| for a model problem of a two-mode atom- 
molecule condensate, and in Ref. jl8| for a description of 
quantum effects on curve crossing. This method is sim- 
ilar to the parametric approximation used in quantum 
optics (see Ref. [[l9)), and in the description of interac- 
tion of electromagnetic and matter waves in Refs. ^0|. It 
is also closely related to the linearization approximation 
in quantum soliton theory (see Ref. [^l] and references 
therein) and to the theory of quantum integrals of motion 
(see Refs. p^] and references therein). The fluctuations 
arc treated here by a linearization of the exact quan- 
tum equations of motion for the field operators. In the 
case of a stable BEC, the resulting linear equations de- 
scribe the motion of quantum oscillators, corresponding 
to the Bogolubov quasiparticles. However, in the case of 
unstable condensates, these equations describe inverted 
oscillators, leading to a rapid growth of the fluctuations. 
This growth remains super-exponential as long as the 
fluctuations do not reach a macroscopic occupation. This 
dynamics offers a more precise picture compared to the 
exponential growth obtained from a Liapunov analysis 
of the classical mean field. Initially, the dynamics can 
be described as a quantum process of spontaneous decay 
of the condensate into the fluctuation vacuum (see Refs. 
p7| , ^8| ), accompanied by the formation of entangled pairs 
(see also Refs. [p3|-p5[ ) . As in the case of fluctuations in 
a stable condensate (see Ref. ^6|), the pairs are formed 
in two-mode squeezed states. 

The formally calculated energy of a Bogolubov quasi- 
particle becomes complex also when the energy of the 
excitation is lower than the energy of a particle in the 
condensate state. This case can take place if the con- 
densate itself occupies an excited translational state, as 
realized recently in the NIST experiments on four-wave 
mixing of atomic condensate waves p^J . The loss of con- 
densate atoms in such a situation was studied in Ref. [^8| , 
using the method of complex scattering lengths, with- 
out taking into account second-quantized properties of 
the atomic scattering. The results of the present theory, 
which take into account these properties, are in agree- 
ment with Ref. [^8| at short collision times [see Eq. ( |6^ ) 
below]. However, at long collision times the effects of 
Bose enhancement lead to a substantial increase of the 
losses [see Eq. (^3|) below]. Notice that the necessity of 
accounting for the quantum effects was mentioned in Ref. 
p8[ itself. These effects should appear under the condi- 
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tions of the NIST experiments [£7| . The present method 
can describe only qualitatively this case of large losses, 
comparable in magnitude to the initial occupation. This 
case may be treated by the incorporation of quantum ef- 
fects into the slowly- varying-envelope approximation de- 
veloped in Ref. [£8]|2!|. 

Quantum effects on the loss incurred by collisions of 
two condensate wavepackets were studied in Ref. [ j30| , 
using an approach similar to the present one. However, 
the results of the present work disagree with those of 
Ref. Q. The losses quoted in Ref. |3(J, which have 
been evaluated only for short collision times, differ also 
from the results of Ref. pq |. These disagreements arc 
apparently due to some u njust ified neglections in Ref. 
|p0| (see discussion in Sec. [IIP ). 

As in the case of fluctuations in a ground-state con- 
densate, the scattered particles are formed here in entan- 
gled (squeezed) states. Therefore the processes consid- 
ered here may open up possibilities for the formation of 
an entangled gas remained after the decay of a conden- 
sate, or after the collision of wavepackets. The entangled 
gas can have many application (see Ref. |26[)^_ 

This paper is organized as follows. Section O considers 
the dynamics of instabilities in a ground-state condensate 
with a negative scattering length. The equations of mo- 
tion for the condensate and the fluctuations are derived 



in Sec. [I A by using a projection onto a coherent s tate . 
The equations for the fluctuations are solved in Sec. II B . 
In the dependence on the excitation energy, the fluctu- 
ations are divided into two classes, of stable and unsta- 
ble modes, which are respectivel y ass ociated with regular 
and inverted oscillators in Sec. II C . The dynamics and 
stati stica l propertie s of t he instabilities are considered in 
Sec. II D , and Sec. HE describes the dynamics of con- 
densate depletion. 



The collision of two matter waves is studied in Sec. Ill 



General equations f or the condensa te and the fluctuations 
are derived in Sec. Ill A. Section IIIB considers elastic 
losses in collisions of two matter waves. The collisions 
of rectangular and Gau ssian wavepackets are treated in 
detail in Sees. Ill C and III D , respectively. 

A system of units in which Planck's constant is ft = 1 
is used below. 



II. GROUND-STATE CONDENSATE 

A. Separation the condensate and fluctuations 

Consider a homogeneous gas of bosons of mass m, in 
a finite quantization box of volume V, described by the 
following Hamiltonian in the momentum representation 
(see Refs. [||), 



2V 



]T U (q) ft (p', t) $ (p - q, t) $ (p' + q, t)] . (1) 



Here ip (p, t) is the annihilation operator for bosons in 
the momentum p-mode, e (p) = p 2 / (2m) is the kinetic 
energy, and U (q) describes the elastic scattering. In the 
model case of a zero-range contact potential, U is a con- 
stant, and can be expressed in terms of the elastic scat- 
tering length a a as 



U = 47ra a /m. 



(2) 



(A modification of the contact-potential model was re- 
cently proposed in Ref. (3lj which allows to avoid ultra- 
violet divergences). The Hamiltonian ([!]) leads to the 
following equations of motion for the annihilation opera- 
tors: 

i$ (p,t) = e (p) ip (p, i) 

u (?) ^ (p'> *) $ (p - q. *) 4> (p' + q. *) ( 3 ) 



q,p' 



Following a general method [^10| let us separate the 
condensate part by a projection onto a coherent state 
1$) defined so that 



4>( P) t)\$)=5 p0 $(t)\$) 



(<D|$) = 1, 



(4) 



where the mean field $ (t) = (-0(p,i)) = (&\ip(p,t) |<&) 
describes the ground-state condensate occupying the mo- 
mentum state p = 0. The annihilation operator can be 
written as 

j> (p, t) = [$ (i) |$)(*| + 1 (P, t)] exp (-iU (0) Nt/V) , 

(5) 

where $ (t) = # (t) exp (iU (0) Nt/V), the operators 
£ (p, t) represent the non-condensate part (fluctuations) 
having the zero expectation value, 



and 



(£ (p,*)) = o > 



(6) 



l*(*)l a + E^ (p>*k a (p>*)> 



(7) 



is the conserved total number of bosons. The exponential 
factor in Eq. (||) shifts the zero energy to the energy of 
non-excited condensate. 

Substitution of Eq. (||) into Eq. (||), combined with 
the property (|Tj|), produces equations of motion for $ (t) 
and £ (p, t). Assuming a small occupation of the fluctua- 
tions, let us keep in these equations only terms of lowest 
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non- vanishing order in £ (p, t) and ^ (p, t), obtaining the 
following equations: 



p 

+^£tf(p)<f(p>*)l(p> *)>*(*) ( g ) 



i£ (p, i) = d (p) | (p, t) + 5 (p) f (-p, t) , 



(9) 



where 



d( P ) = e(p) + U(p)n , g(p) = U(p)<£ 2 /V, (10) 

and no = |$| 2 /V is the condensate density. This approx- 
imate method, that leads to linear operator equations for 
the fluctuations, is similar to the parametric approxima- 
tion used in quantum optics (see Ref. |l9]| ). It has been 
applied to the study of a molecular BEC dissociation in 
Refs. |l^,|l^]), and is valid as long as the depletion of the 
condensate is negligibly small. An a pplica bility criterion 
for this method is presented in Sec. HE below [see Eq. 

The present approximation is gapless [see Eq. flic 
below] and, therefore, non-conserving (see Ref. B), al- 
though the average number of the atoms is conserved 
[see Eq. (0) above] . 



B. Stable and unstable modes 

Equation (JoJ) combined with the corresponding equa- 
tion for £t (— p, t) provide a set of two linear operator 
equations (similar to the set considered in Refs. |l7| , |l8| ] ) . 
This set can be solved exactly, using common methods 
reducing it to an eigenproblem, by assuming the time- 
independence of d (p) and g (p) . This assumption is cor- 
roborated by Eq. (^), which shows that the variation of 
$ is of the second order in £ (p, t) and therefore can be 
neglected. Indeed, the validity of the parametric approx- 
imation requires a small variation of the condensate den- 
sity, and the dominant part of the condensate phase is 
accumulated in the exponential factor in Eq. (g). 

As in the case considered in Ref. [[lTj , the form of the 
solution depends on the ratio of the detuning d (p) to the 
coupling strength g (p) . All momentum modes observing 



\g (p)\>d (p) , or e (p) + 2U (p) n < 



(11) 



are unstable, and the corresponding solutions attain the 
aperiodic form 



e(p,*) 



where 



cosh (A (p) t) - i^EL sinh ( A (P) *) S (P> °) 
- ? |Msinh(A(p)tKt(-p,0), (12) 



A(p)= V\g(p)\ 2 -dH P ) 

= y/e(p) (-2U(p) no-e(p)). 



(13) 



The condition ( |ll[ ) can be obeyed only when U (p) < 0, i. 
e., with a negative scattering length, and for excitations 
with an energy e (p) < 2\U (p) \uq. Otherwise, whenever 
\g (p) ] < d(p), the modes are stable and the solutions are 
oscillatory functions of the form 



£(p,*) 



cos(|A(p)|i) 
, g(p) 



l 4^sin(|A(p)|t) 



|A(P) 



A(P) 

Bin(|A(p)|t)f (_ p ,o). 



£(p,0) 
(14) 



In a case of intervening mode (A(p) = 0) the solution 
is a linear function of time 

i (p, t) = [1 - id(p)f] £ (p, 0) - %(p)tf (_p, o) . (15) 

This case can take place for the spurious (Goldstowne) 
mode with a vanishing excitation energy (e(p) = 0, see 
Refs. |fl), or at e(p) = 2\U(p)\n . 



C. Regular and inverted oscillators 

The stable modes can be treated by the Bogolubov 
transformation. It is therefore instructive to establish a 
connection between the present and the Bogolubov meth- 
ods. The equations of motion (^|) can be derived from the 
following Hamiltonian constituent involving two counter- 
propagating modes, 

h ( P ) = d (p) [|t ( P , t ) i (p, t) + |t (_ Pj t) I (-p, t) 

+5* (P) f (p, *) £ (-P, f) + .9 (p) f (p, *) I 1 (-P, *) ■ (16) 
The transformation 

li,2 (p,t) =2-Va (£(p,t)±£(-p,i)) exp (-?- axg g(p) 

(17) 

introduces a new pair of field operators, £i (p, t) and 
£2 (p, t) , whose equations of motion are decoupled by sep- 
arating the Hamiltonian into H (p) = H\ (p) + if 2 (p), 
where 



d£|6 - i (-i) j | 5 | (ell] + 64) > i = 1. 2 - (i 8 ) 



(Hereafter in this subsection the arguments p and t are 
dropped out, for simplicity's sake.) A system described 
by the Hamiltonian ( |18|) was solved in Ref. |l7j , and the 
solutions have a form similar to Eq. (|l2|), or Eq. (H), 
according to the ratio of d to 5. 

Let us introduce canonically conjugate operators for 
the coordinates Qj and the momenta Pj, such that £j = 

(<5j + i-Pj) , writing out the Hamiltonian ( |l8| ) as 
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d+(-iY\g\ 



p 



d-(-iY\g\, 



d 

2' 



(19) 



If \g\ < d, this Hamiltonian describes a harmonic oscilla- 
tor and can be transformed into the standard form repre- 
senting free Bogolubov quasiparticles [the total transfor- 
mation, incorporating (j!7|) , being indeed the Bogolubov 
transformation] . The motion of the harmonic oscillator 
is periodic, in agreement with the periodic behavior of 
the solution (pj). 

A different situation takes place if the condition Jll| ) 
holds. In this case one of the coefficients (preceding P? or 

Qj) is negative, and the Hamiltonian (|l9| ) represents an 
inverted oscillator. It cannot be transformed to a form 
representing free quasiparticles, and therefore a Bogol- 
ubov transformation cannot exist, in agreement with Rcf. 
H . The infinite aperiodic motion of the inverted oscilla- 
tor is exhibited in Eq. (|l^) and in the instability of the 
mode. 

In a case of intervening mode (\(p) = 0) Eq. (Jl^) gives 
d{p) — \g{p)\, one of the coefficients (preceding Q? or P?) 
vanishes, and the Hamiltonian ( JTgj ) represents a motion 
without restoring forces or an oscillator with infinitely 
large mass, respectively (see Ref. ||). Such systems 
demonstrate linear time-dependence of Pi or Qj, respec- 
tively, in agreement with the solution (Uq). The Bogol- 
ubov transformation does not exist in this case as well. In 
the case of the Goldstowne mode (p = 0) the Hamiltonian 

(16) has already the form (|l8|), and the transformation 

(17) is not necessary. This mode brings only one term 
Q 2 (0, i) into the total Hamiltonian, in agreement with 
Rcf. |3^] (up to the exchange of the coordinate and the 
momentum). 



D. Dynamics of instabilities 

Consider the dynamics of fluctuations in an unstable 
BEC, starting from the state of a pure condensate 

(£ f (P, o) i (p, 0)) = (I (p, 0) £ (-p, 0)) = 0. (20) 

This problem is not only of a speculative interest; it actu- 
ally describes the recent JILA experiment Q], in which 
the elastic scattering length was varied very fast from the 
zero value (of the ideal condensate) to a negative value 
(of an unstable condensate) . Given the initial conditions 
(pp|), Eq. ( |l2| ) leads to the following dynamics for the 
occupation of the unstable modes, 

(f (p, t) i (p, t)) = sinh 2 (A (p) t) , (21) 

and for the two-particle correlations, 
, 9 (p) 



(£(p,t)£(-p,t)> 



A(p) 

d (p) g (p) 



A 2 (p) 



cosh (A (p) t) sinh (A (p) t) 
sinh 2 (A (p) t) . (22) 



It is not surprising that this dynamics is similar to 
the dynamics of spontaneous dissociation of a molecular 
condensate, studied in Ref. fTi| , as the non-condensate 
modes can be described in terms of the fields £1,2 (p, t) 
[see Eq. jl^)] by the Hamiltonian ( [l8|) , which is equiv- 
alent to the one considered in Ref. JL7||. By analogy to 
the results of Ref. [O], the non-condensate atoms are 
formed in squeezed states involving atom pairs of oppo- 
site momenta. These states are similar to the two-mode 
squeezed states of fluctuations in a stable BEC studied 
in Ref. ^(|. The states are perfectly squeezed and have 
zero variance in the relative number of particles with mo- 
menta p and — p, being optimally entangled in the en- 
tropic sense. The uncertainty of the generalized coordi- 
nates 



Xj (p, t) = 2- 1 ' 2 [fc (p, t) exp (%9j) + £} (P, t) exp (-ity) 

(23) 

where the operators £j (p,t) (with j — 1,2) are defined 
by Eq. (|l^), attains the maximal and minimal values 



(A >2 (p,i)) ± = 



d 2 (r>\ 

cosh 2 (A (p) t) + — f{ sinh 2 (A (p) t) 



A 2 (p) 



1/2 



±||| S inh(A(p)t) 



At»l 



\g(p)\ r 2\t 

A» r -2Xt [M) 



at the two orthogonal values of the respective phase an- 
gles 9 j, 



J± =±(-1) J 



1 



arctan 



dip) 
Hp) 



tanh (A (p) t) 



(25) 



A common method for studying BEC stability is the 
Liapunov analysis (see e. g. Ref. 0). This method pro- 
duces equations for small perturbations Q (p, t) to the 
mean field. These equations are equivalent to Eq. (|^), 
in which the operators £ (p, t) are replaced by c-number 
functions £(p,i). In the case studied here one obtains 
for the Liapunov exponents the expression A given by 
Eq. ([l3"|). This type of analysis attributes an exponential 
gain to the small perturbation, |£(p,£)| 2 ~ exp(2A£), 
when the condition (11) is obeyed. This result is in qual- 
itative agreement with the quantum analysis used here, 
which attributes an even faster (super-exponential) gain 
to the fluctuations [see Eq. (21)], as 



dt 



In 



(<f (P, t) i (p, t)>) - 2A (p) coth A (p) t > 2A (p) . 

(26) 



The dynamics of the fluctuations has been studied in 
Rcf. |17| . It is shown there that the quantum (squeezed) 
state always grows faster than a classical (coherent) 
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state. Therefore the Liapunov analysis concerning classi- 
cal fields underestimates the fluctuation growth, becom- 
ing correct only when At ^> 1, where the fluctuations are 
large enough to be considered as classical ones. 

A formal application of the Bogolubov transformation 
at U < (although incorrect, as had been observed by 
Bogolubov himself; see Ref. ||) produces quasiparticles 
with imaginary energies iA, which correspond to expo- 
nential growth. However, the transformed field operators 
do not obey the proper Bose commutation rules. 

The homogeneous boson gas with a negative scatter- 
ing length considered above always has unstable modes 
due to the continuous nature of the excitation spectrum. 
The situation changes in the case of a trapped gas having 
a discrete excitation spectrum. Following Ref. [Q, con- 
sider A (p) as a function of the condensate density no for 
a fixed value of p = p m chosen so that e (p m ) = £m, where 
e m is the energy of an excited state of the trapped con- 
densate. At low condensate density A (p m ) is imaginary, 
corresponding to the well-known stability of a trapped 
condensate with a negative scattering length at low oc- 
cupations (see Ref. Q). At n > |e m /|J7(po) I the value 
of A (pm) becomes real, and therefore describes an in- 
stability of the condensate (leading to its collapse). Al- 
though the present approach does not take into account 
the inhomogeneity of the condensate and the variation 
of its size, the dependence of A (p m ) on no [see Eq. (f 3)] 
is in qualitative agreement with the one obtained from 
the numerical calculations of Ref. |33| for the monopole 
mode. 



section. The statistical weight of intervening mode with 
A(p) = is negligible compared to the continuums of 
the stable and unstable modes. Unlike a case of sta- 
ble condensate, where the unrestricted behavior of the 
Goldstowne mode requires a special treatment (see Ref. 
|p2| ), in the present case the effect of linearly-growing 
intervening modes is negligible compared to the effect of 
exponentially-growing unstable modes. 

Substitution of p — p s cos | allows us to express / U nst 
in terms of the Weber functions E v (z) (see Ref. p3): 



7T 

f 9 

lunst— J cos — sinh (r sin 9) d6 



= «2 [Ei/2 (it)+E_ 1/2 (it)] , 



(30) 



where r = 2£/£nl, and £nl = (l^l^o) is the nonlinear 
interaction time (see Ref. [^9] ) . Approximate expressions 
for I uns t, holding at small and large values of t, have the 
following respective forms, 



r^o 2 v / 2r(7/4)r(5/4) 



T ~ 1 /Lt 3 



(31) 



t^oo 2 V t 16r 



The contribution of the stable modes can be evaluated 
by substituting p = p s cosh | , as a result of which 



E. Dynamics of the condensate 

Let us consider now the dynamics of condensate deple- 
tion. In the limit V -> oo, Eqs. (§!]), (|| ), and (§) lead 
to the following equation for the condensate density, 



no 



rnp s U Uq 
2^2 



where 



'stab 



mp s 



mp s 



(Ainst ~\~ -^stab) ; 

a sinh (2A (p) t) , 

P Hp) dp 

2 sin(2|A(p)|i) 



(27) 



(28) 



|AO) | 



dp 



are the contributions of the unstable (p < p s ) and stable 
(p > p s ) modes, respectively, 



p a = y/-im\U\n 



(29) 



is the boundary between the stable and unstable modes, 
and the contact potential (0) is used hereafter in this 



/ 



-fstab= / cosh — sin (r sinh 9) d9 



= l\l^ £T " l \ t El/2 {lT) + E - 1/2 (iT) ] ' (32) 
Its limiting values can be approximately expressed as 
7T ( I It t 



-^stab 



r^o 2 \y/¥f Mir V2r(7/4)T(5/4) 



(33) 



stab 



t^oo 16r 



Finally, the rate of condensate depletion can be expressed 
as 

mp s U 2 riQ T 
n = : 



Att^/ttt 

-4W — a^ngexp (87r|a a |n t/ra) , 



(34) 



and the loss of condensate density attains the limiting 
forms 
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An « 4noy/\a a \ 3 n T = ?,J^-a 2 a nl, t < t NL (35a) 



An « 2e r n Vl a a| 3 ™o/ 7 



a a n exp (87r|a a |noi/m) , t 3> tNL (35b) 



at the two distinguishable loss regimes. At the regime 
of non-stationarity t -C tNL the depletion is mostly due 
to an occupation of the stable modes, which carry much 
more statistical weight than the unstable modes [I s t a b 3> 
^unst, see Eqs. (|lj) and (33)]. The statistical weight is 
infinitely large, leading to a singularity in the depletion 
rate ([Mj), at t = 0. This singularity is absent in Eq. 
(|35|), and can be eliminated by using of an interatomic 
interaction with a finite range in place of the contact 
potential (||) . At t > the occupation of the modes with 
p > y/m/t begins to decrease due to their oscillating 
behavior, leading to decrease of the depletion rate. At 
the regime of Bose-Enhancement t ^> <nl the depletion 
is dominated by the exponentially increasing occupation 
of the unstable modes [7 s tab <C iunstjSee Eqs. ( |3l"| ) and 
(|33|)1. The approximate expressions (35) are compared in 
Fig. [| to results of numerical integration of Eq. (34). 

One should be reminded that the parametric approxi- 
mation, on which all the analysis above rests, is applica- 
ble only as long as Ano <C tiq. In the case of one excited 
mode JlTj , this approximation is in good agreement with 
numerical results while Ano/no < 0.2. In the case of a 
dilute condensate (|a a | 3 no -C 1, sec Rcf. Q), the validity 
range 



*/*nl < -|ln(|a a | 3 n ) | 



(36) 



includes a region t > tNL where the exponential Bose- 
enhancement factor can reach a rather large value, of the 

order of (|a a | 3 7to) ^ 2 . 

Equation ( (35a| ) predicts small depletion at times t < 
tNL (with An /n < 4y/\a a \ 3 no <C 1, given the dilution 
parameter |a a | 3 no is small). Therefore, the development 
of instabilities in BEC, including collapse, starts with 
a delay time i ~ tNL after turning on the interatomic 
interactions, although the collapse itself is a nonlinear 
phenomenon that cannot be treated by the parametric 
approximation. Under the conditions of the 85 Rb exper- 
iment 14 (where no w 10 12 cm -3 and a a ~ — 20nm), the 
expected delay is tNL ~ 5ms (cf. Fig. |l|), in agreement 
with the experimentally observed delay time for the col- 
lapse. In this case the parametric approximation is valid 
as long as t < 3tNL, allowing for values of the exponen- 
tial Bose-enhancement factor in Eq. ( |35b| ) exceeding 100. 
It should be noted that the condensate dynamics can be 
sensitive to initial correlations and the spatial inhomo- 
geneity neglected here. 



III. COLLISIONS OF TWO CONDENSATES 



A. Separation of condensates and fluctuations 

The approach developed in Sec. || above can also be 
applied to a collision of two free counterpropagating BEC 
waves with momenta ±po (see Ref. |2^|). In this case the 
annihilation operator should be represented as 



V^(p,t) 



4>+ (t) |$+)($+ 1 exp (-2zS_ (t)) 
+$_(t) !$_)($_ | exp (-2*S + (t)) 
+i(p,t)exp(-iS + (t)-iS- (t)) 



(37) 



where 



2V 



J \$±{t>)\*dt> +^Nt+\e{p )t. 



(38) 

Here <&± (t) describe the two condensate fields with mo- 
menta ±p , respectively, the coherent states |$±) are 
defined by 



^ (p, t) |$±) = *± (t) exp (-2iS T (t)) <5 ±P0 p|$±) 

(<D±|$±) = 1, 
and the conserved total number of atoms is 



(39) 



N = i$+ (*) i 2 + (*) i 2 + £<f (p» *) z (p» *)) • ( 4 °) 

p 

Substitution of Eq. (|37j) into Eq. (|J) gives the following 
equations of motion for the fluctuations, 

i£ (p, t) = d 2 (p) i (p, t) + 52 (p) f (-p, t) + D, (41) 
where 



d>2 (p) =e(p) -e(po) + 



^(|p-Po|)-if/(2p ) 



n i 



C/(|p + Po|)-2^(2Po) 



(42) 



g 2 (p)=[U(\p-p \) + U(\p + p \)}$ + $-/V, (43) 

and n± = | < i ) ±| 2 /V are the densities of the condensate 
waves. The term 

D = {$*_ (t) $+ (t) [C/ (2 Po ) + 17 (|p - po|)] f (p - 2 Po , t) 
(t)f/(|p-p |)e t (-p + 2p ,t)} 
xexp{2i[5+(t) - S_(t)]} 
+{$; (t) 4>_ (t) [[/ (2p ) + C/ (|p + po|)] i (p + 2p , t) 

+$ 2 _(t)(7(|p + po|K t (-p-2po,*)} 
xexp{2i[5_(t) -S+(t)]} (44) 



G 



describes couplings of fluctuations with momenta differ- 
ent by 2po. Whenever 



B. Elastic-scattering loss in collisions of two matter 
waves 



pi > 4w£/ v /rEj 



(45) 



the disparity of the energies of such fluctuations far ex- 
ceeds the coupling strengths, and the term D can be 
neglected in Eq. (f4l|). This allows us also to neglect cor- 
relations of such fluctuations in the equations of motion 
for the condensate fields, which then attain the form 

**± (*) = p £ u dp - poi) d (p> f ) £ (-P- (*) 



~ 2 c dp - poi) <e f (p, t) I (p, m± (t) (46) 



Equation (|4l|) with D neglected has solutions of the form 
Eqs. (|I^), (|14j), or (|l5|), according to the ratio of |<7 2 | to 
d,2i in which <?, d, and A are respectively replaced by <? 2j 
d 2 , and 



A 2 (p) = ^152 (P) I 2 



:(p). 



(47) 



The condition of instability <i 2 (p) < |<?2 (p) | can be 
obeyed now even for positive values of the scattering 
length (and hence of U), as d 2 contains a negative term 
— e(po) [see Eq. (|4^)]. In the case of a contact potential 
(0) an instability takes place when 



\e{p)-e(p ) + ~U(n++n-)\ <2Uy/n^nZ, (48) 
or when 



where 



P± 



P. < P < P+, 



m/7 (n_)_ + n_) ± 4m[/^/n+n_ 



1/2 



(49) 



(50) 



Whenever condition ( f45| ) is obeyed, the momenta p± are 
given approximately by 



mil n + + 7i- 2mU 

P± ~ Po ^ ± y/n+n. 

Po 2 po 



(51) 



The instability interval can include the region p > p a . It 
may appear that this fact can lead to the absurd con- 
clusion that a ground-state condensate (po — 0) with 
a positive scattering length is unstable. However, the 
analysis above is valid only when the condition (B^) is 
obeyed. Otherwise, the terms neglected in Eqs. (|4l|) and 
( [itf ) have to be taken into account, leading to stability 
of the condensate at po — 0. 

As in the case of a ground-state condensate, the scat- 
tered atoms with momenta p and — p are correlated, 
forming optimally entangled and optimally squeezed 
states [see discussion following Eq. (|22|)]. 



Consider now the elastic-scattering loss of atoms from 
the two colliding counterpropagating waves bounded only 
by the normalization box. Substitution of Eqs. ( fH| ) and 
( p2| ) (with the replacement of d, g, and A by <i 2 , g2, and 
A2 , respectively) into Eq. ( fl6| ) leads to the following equa- 
tions for the two condensate densities, 



n± 



mpoU 2 n + n- 



I-) 



(52) 



where 



■Minst — 



mp 



2 sinh(2A 2 (p)t) j 

p \ — dp 



A 2 (p) 



7_ = 



mpo 



I + = 



mpo 



2 sin(2|A 2 (p)|t) 
P IA2 (P) I * 

p2 sm(2\\ 2 (p)\t) 



(53) 



|A 2 (P) 



dp 



represent the contributions of the unstable (J un st) an d 
stable (I + + J-) modes, respectively. 
The substitution of 



P 



mU (n + + n_ ) + imU y/n + ri- cos ( 



1 1/2 



allows us to express I U nst in terms of the modified Struve 
function L v (z) (see Ref. |34j]), as 

iumrt = / M— sinh (r sin 9) » ttL (t) , (54) 
J Po 



where r = 2t/t^i n and here £nl = (2C7 ^/n + n_) . Ap- 
proximate expressions for 7 uns t at small and large values 
of t have the respective forms 



Ainst ~ 2t, /, 

T— >0 



unst ~ A/ ^ 



(55) 



Given the condition tp\jm <C 1, the substitution of 
IA2 (p) I ~ p 2 / (2m) produces the estimate 



1 

po 



(56) 



The integral 7_ can be estimated under the same condi- 
tion by retaining only the first term in the Taylor series 
for sin (2|A 2 (p) \t), giving J_ « 2pgt/ (3 m ). 

Given the opposite condition, tp\jm 3> 1, the con- 
tributions of the stable modes I± can be evaluated by 
substituting 



7 



p = [pq — mil (n+ + n_ ) ± 4mUy/n+n- cosh ( 
resulting in 



I 1/2 



/+ = 



d6^^sin(Tsinh0). 



Po 



(57) 



where 9 + = oo and 

ln[pg/(2mCV™+^)])- 
stable modes, given by 



p{6-) = (with 6- w 
The total contribution of the 



/++/_ W7r(/ (r)-io (r)) : 



(58) 



where Jo (z) is the modified Bcssel function (see Ref. 
p4j), can be approximately expressed as 



7T - 2t, 
2/t, 



ra/p 2 , < i < ^nl 
t Xnl- 



(59a) 
(59b) 



Combining Eqs. (|52|), (|54|), and (58), one obtains the total 
rate of the condensate loss 



n± 



mpoU 2 n + ri- 



Io(r) 



(60) 



Unlike the case of ground-state condensate, character- 
ized by a single time-scale tjMLj an additional time-scale, 
tti/pq, appears in the present case, leading to appearance 
of a third regime of the condensate loss. At short times, 
t <C £nlj the loss is dominated by the occupation of the 
stable modes due to their higher statistical weight. 

At the regime of non-stationarity tp\jm <C 1, the 
wavefunction of the scattered particles does not have 
enough time to become stationary, and the estimate (|5^ ) 
gives the loss rate 



n± 



m!7 2 n+n_ / m 
7T V 2~7rf' 



(61) 



The singularity at t = has the same nature as in the 
case of the ground-state condensate with a negative scat- 
tering length [see discussion following Eq. (p5|)] . 

At the regime of stationarity, ra/p 2 , «(« £nl, the 
loss rate becomes time independent, and given by 



n± 



ic 2P0 



(62) 



in full agreement with theory p8[ . 

At the regime of Bose-enhancement, i ^> £nl, the loss 
is principally determined by the occupation of unstable 
modes which increase exponentially due to the effect of 
Bose enhancement. The loss rate can then be approxi- 
mately expressed as 



n± 



(2a a f 2 p (n+n_) 3 / 4 



'mt 



exp(2t/t 



nlJ 



(63) 



The approximate expressions (|6l]), (|62|), and (^) are 
compared in Fig. |^ to results of numerical integration in 
Eqs. @. 



The contribution of stable mo des decreases at long 
times as t increases [see Eq. ( 59b )]. In order to explain 
the decrease, let us consider the asymptote of I + + 7_ 
using, in the spirit of scattering theory [pa] , the substi- 
tution 



hm sin (2jA 2 (p) \t) /|A 2 (p) | = 7r<5 (A 2 (p)) 
<*(p-P+) tf(p-p_) 



^|A 2 (P+)| £|A 2 (p_) 



(64) 



This leads to 



2mp 



i\X 2 (p + )\ £|A a (p_)| 



(65) 



In the absence of many-body effects ([/ = 0) , one ob- 
tains |Aa (p) | = |p 2 — p 2 |/(2m), expressing the regu- 
lar dispersion law for free particles. Also, p± = po, 

^|A 2 (p±) I = Po/m, and hence + J_ = 7r, which 

also hold for small t [see Eq. ( |59a| )]. Therefore the scat- 
tering process can be described by a stationary scattering 
theory at any time t 3> mjp\ (see Ref. |3q[). Contrarily, 
in the case of a finite t7, many-body effects modify the 



dispersion law, so that 



U + I- 



0. 



^|A 2 (P±) 



00, and hence 



C. Collisions of rectangular wavepackets 

The discussion above is related to the scattering of 
two condensate waves bounded only by a normalization 
box. A more realistic yet simple model is that of a half- 
collision of rectangular wavepackets. Consider a homo- 
geneous condensate occupying a region of length b along 
the x axis (leaving the results independent of size and 
shape in the yz-plane). At the moment t = a Bragg 
scattering pulse is applied, creating a wavepacket of den- 
sity n+ and momentum 2po directed along the x axis, 
while the rest of the atoms, with density n_, remains in 
the p = wavepacket. Losses take place only at the over- 
lapping parts of the two wavepackets. The overlapping 
time tovor varies from to t co \, between the leading and 
trailing edges of the moving wavepacket, where 



tr 



mb 
2p^ 



(66) 



is the total duration of the collision of the two wavepack- 
ets. The average density loss can be expressed as 



tool tovcr 

(An±) = J dt OVCI J dth± (t) 





where n± (t) is defined by Eq. 



(67) 
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At the regime of non-stationarity, t co pQ/m < 1, or 
bp <C 1, the loss rate is given by Eq. (pi), and 



(An±) 



.32 



-a a n+n_ 



(68) 



At the regime of stationarity, m/p 2 <C i C ol -C £nl, the 
loss rate is determined by Eq. (62), and the loss has a 
form 



(An±) ps 47ra 2 &n_|_n_ 



(69) 



The same expression for the loss can also be obtained by 
using the method of complex scattering length ps| ]. 

At the regime of Bose-enhancement, t co \ ^> £nl, the 
loss rate Eq. (63) gives the average density loss 



(An±) 



5/2 

Po 



32ir 2 a a /2 b 3 / 2 (n+n_)V4 



exp (2tcoi/t NL ) (70) 



The approximate expressions (|6q), (69), and ( pXjj ) are 
compared in Fig. to results of numerical integration 
in Eqs. © and ©. 

The parametric approximation and, hence all the anal- 
ysis above, are applicable only as long as (An±) <C 
. Using Eq. (f70|), this condition can be rewrit- 
ten in the form of a limit on the duration of collision 
consistent with the approximation: 



1 327r 2 ayV/ 2 (n + n_) 3 / 4 

^col < ^^NLhl gy^ . 

Po 



(71) 



For example, in the NIST Na experiment (where 
?i + rj n_ 10 14 cm~ 3 , a a ps 3nm, b ps 20/im, and 
Po ~ 10 5 cm _1 ; see Ref. ^9|), the parametric approx- 
imation is applicable while t co \ < 3£nl, which allows 
for values of the exponential Bose-enhancement factor 
in Eq. ( [70| ) of more than a 100. This validity range sub- 
stantially exceeds the range of applicability of the com- 
plex scattering- length method, t co \ <C £nl- However, in 
the four- wave- mixing NIST experiment p7[ , t co \ ps 7£nl, 
which requires taking into account a considerable deple- 
tion of the condensates during the collision. This task 
may be performed by combining the present approach 
with that of the slowly-varying-envelope approximation 
ppl. 



ii(r,t) = ^(r,t)+N J dV(7(r-r') [|^(r',t) | 2 |(r,i) 
+<p* (r', t) if (r, t) I (r', t) + <p (r, t) <p (r', t) f (r', i)] , (72) 

where </? (r, t) is a unit-normalized wavefunction describ- 
ing particles in the colliding wavepackets, N is the total 
number of particles in the wavepackets, and U (r) is the 
interatomic interaction potential. 

The first-order perturbative solution of Eq. ( frgh can 
be expressed in the form 



& (r, t) = j d 3 r' J dt> [A (r, t; r', t') £ (r', 
o 

+B(r,t;r',t')f (r',t')]. (73) 
Here the zero-order solution 

lo(r,t)=^C(p ) 0) X p(r ) t) (74) 
p 

is expressed in terms of free (plane) waves 



Y P (r,<) = V- 1 / 2 exp(zpr- l |-<) (75) 



and the annihilation operators £ (p, t) for the fluctuation 
in the momentum representation. As the following equa- 
tions turn out to be independent of the first term in Eq. 
([73|), including the function A (r, £; r', t') , it is not neces- 
sary to specify that function. The second term includes 
the function 

B (r, t; r', t') = -*7V f d 3 r" ]T Xp (r, t) X * p (r",t') 
J p 

xU{r" -r')ip(r",t')<p(r',t'). (76) 

Considering scattering into the vacuum of fluctuations, 

such that (f (r,i)f (r',t')) = (f„ (r,t)|„ (r',i')) = 0, 
and therefore 

(| (r, t) |o (r', f )) = £ X P ( r < *) *P ( r '' • ( 7? ) 



D. Collisions of Gaussian wavepackets 

The consideration of wavepackets of a more realistic 
shape requires the solution of Eq. (^d]) for time-dependent 
o?2 and <?2- This solution can be obtained by using a per- 
turbation theory, which is applicable as long as t <C £nl, 
and does not describe the effects of Bose enhancement. 
The evaluation is made simpler in the coordinate repre- 
sentation, in which case the linearized equations of mo- 
tion for the field operators of the fluctuations £ (r, t) have 
the form 



one can express the number of the lost atoms in the form 

-J 



AN = -J d3r (ti(r,°°)Zi (r,oo)) 

OO OO 

= - / d 3 r J d 3 r 1 J d 3 r 2 J dt x J dt 2 B* (r,ao;ri,ii) 

00 

xB (r,oo;r 2 ,i 2 )£x P (ri,ii)x* (r 2 ,t 2 ) 
p 

(78) 



2Af 2 ImT ( o ) 
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where 



oo oo 



T ( 2) = J dt" J dt' J d 3 r[d 3 r' 2 d 3 r'{d 3 r%<p* (r'/, t") 
o o 

x^* (r'i, t") U {v'{ - r'i) G (r'{, r%, t"; r^, t>) 
xU(r\-r' 2 )^(v\,t') V (r' 2 ,t'). (79) 



Here 



G(rlr>>X;r[y 2 ,t>) = -i £ X Pl (4,t")x P2 (4X) 



XX* P1 (r[,t')x* m (r 2 ,t')#(t"-t') 



(80) 



and -d (t) is the Heavyside step function. As one can see 

G(r", r 2 ,t"; r[, r' 2 ,t') is a time-dependent Green function 

i 2 \ 

for two free particles. Therefore T 00 is nothing else but 
the second-order term for the diagonal matrix element 
of the transition amplitude, obtained by time-dependent 
perturbation theory for a collision of two particles, each 
in the state ip(r,t). If the size of the wavepacket obeys 
b l/po, the result can be reduced to a form given by 
the stationary scattering theory (see Ref. |35|), which is 
proportional to the imaginary part of the elastic scatter- 
ing length, in agreement with Ref. ]2S(| . Consider two 
colliding Gaussian wavepackets described by the wave- 
function 



<p(r,t) 



where 



F+ ( r - ^-t) e lpor + F_ (r+^i ) e~* p " r 



x exp 



2m > ' 



(81) 



F± (r) = n- 3 / 4 ^/N(b x b y b z y 1/2 

and N± is the number of particles in each wavepacket. 
Substituting Eq. (|8l|) into Eq. (|79|) one can obtain 



b x b y b z m 



where 



teal = m 



bl + bl + b\ 



-1/2 



(83) 



(84) 



the one presented here, with solutions of the form of Eq. 
(|l2|) published previously in Ref. [0. However, unjusti- 
fied neglections of substantial terms in the Hamiltonian 
lead to an expression for 172 which is smaller by a factor 
of four than the one given here by Eq. ( [l3"| ) . This neglec- 
tion should therefore result in an underestimation of the 
condensate loss. 

The above-mentioned work, Ref. considers colli- 
sion of Gaussian wavepackets to first-order perturbation 
theory only, which do not describe the effects of Bose 
enhancement. Finite analytical expressions of the kind 
presented here were not derived there, and the results of 
numerical calculations disagree with results of the com- 
plex scattering length method |^8| even after the cor- 
recting for the missing factor four in g 2 are taken into 
account. In contrast, the analytical expression Eq. ( |8^ ) 
of the present paper, based on the same approximations, 
is in an agreement with the complex scattering length 
method. 



IV. CONCLUSIONS 

Fluctuations in BEC are treated by using a paramet- 
ric approximation for the exact quantum equations of 
motion for the field operators. Depending on the ex- 
citation energy, the fluctuations are divided into stable 
modes, which can be represented as Bogolubov quasipar- 
ticlcs (field oscillators), and unstable modes, which can 
be represented as inverted oscillators and demonstrate a 
rapid growth. Initially the dynamics of unstable modes is 
a quantum process of spontaneous decay into the vacuum 
of the fluctuations with a super-exponential rate. 

This approach is applied to two physical problems. 
The first one involves the dynamics of a BEC after a 
fast change of the scattering length from a non-negative 
to a negative one. In this case , the experimentally ob- 
served Q delay of the condensate collapse is described. 
The second problem involves the elastic scattering losses 
on collision of two BEC wavepackets. At short colli- 
sion times [see Eq. (|(32|)1 the losses obtained here are in 
agreement with those of the method of complex scat- 
tering length (28) , with specific reference to the colli- 
sions of non-confined waves, of rectangular and of Gaus- 
sian wavepackets. At longer collision times [see Eq. (joS))] 
quantum effects lead to Bose-enhancement of the elastic 
scattering losses. The atoms lost from the condensates in 
both cases form correlated pairs in optimally entangled 
and optimally squeezed states. 



is the collision duration. Equation fl85| ) coincide with 
the expression which can be obtained by using of the 
complex scattering length method |2^[ and neglecting the 
depletion of the wavepackets. 

The recent preprint p0| | , posted during preparation of 
the present paper, treats the loss of atoms from collid- 
ing BEC wavepackets by using an approach similar to 
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t(ms) 

FIG. 1. Time-dependence of the loss of condensate den- 
sity Ano under conditions of the 85 Rb JILA experiment Jl^] 
(no ~ 10 12 cm~ 3 ). The results of numerical integrat ion (solid 
l ine) are compared to the approximate equations (35a) and 
(35b) (dashed line, left and right segments, respectivelly). 
The dot-dashed line represents the condensate depletion rate 
no given by Eq. 
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b(nm) 




FIG. 2. Time-dependence of the condensate loss rate 
h± under conditions of the Na NIST experiment |27]] 
(n± ~ 10 14 cm -3 ). The results of numerical integration in 
Eqs. ( p3| ) (dot-dashed line) are compared to the approximate 
equations @, (||), and (|3|) (dot-dot-dot-dashed line, left, 
central, and right segments, respectivelly) . The solid line 
represents the average density loss Ano for the collision of 
rectangular wavepackets calculated by using of numerical in- 
tegration in Eqs. ( ^ ) and do7l) , in comparison with the ap- 
proximate expressions (jiiil) , (|6it) , and ( |7o| ) (dashed line, left, 
central, and right segments, respectivelly). 
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